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ABSTRACT 

Using a new large-scale (~ 0.75 Gpc) 3 hydrodynamic cosmological simulation we investigate the 
growth rate of supermassive black holes in the early universe (z 4.75). Remarkably we find a clear 
peak in the typical Eddington ratio (A) at black hole masses of 4 — 8 x 10 7 M Q (typically found in halos 
of <~ 7 x 10 11 — 1 x 1O 12 M0), independent of redshift and indicative that most of BH growth occurs in 
the cold-flow dominated regime. Black hole growth is by and large regulated by the evolution of gas 
density. The typical Eddington ratio at a given mass scales simply as cosmological density (1 + z) 3 and 
the peak is caused by the competition between increased gas density available in more massive hosts, 
and a decrease due to strong AGN feedback that deprives the black hole of sufficient gas to fuel further 
rapid growth in the high mass end. In addition to evolution in the mean Eddington ratio, we show 
that the distribution of A among both mass-selected and luminosity-selected samples is approximately 
log-normal. We combine these findings into a single log-normal fitting formula for the distribution of 
Eddington ratios as a function of (Mbh,z). This formula can be used in analytic and semi-analytic 
models for evolving black hole populations, predicting black hole masses of observed quasars, and, in 
conjunction with the observed distribution of Eddington ratios, can be used to constrain the black 
hole mass function. 

Subject headings: quasars: general — galaxies: active — black hole physics — methods: numerical — 
galaxies: evolution 



1. INTRODUCTION 

It has been well established that supermassive black 
holes are present in the c enter of most galaxies 
(|Kormendv fc Richstonelfr995l ) . and th at they are corre- 
lated with the properties o f their hosts (iMagorrian et all 
19981 iFerrarese fc Merrittl I21M iGebhardt et all 120001: 
Tremaine et all 120021 : iGraham fc Driver! 120071 1 These 
correlations provide strong evidence that the growth of 
a black hole and the evolution of its host galaxy directly 
influence one another, such that black hole growth is a 
important aspect of understanding galactic evolution and 
vice versa. 

In general, the link between black hole and galactic 
evolution is attributed to some form of quasar feedback 
dBurkert fc SilU200lHSazonov et aT1l2004l:lSpringel et all 
2005bt iChurazov et all 120051: iDi Matteo et al l 120051: 
Bower et all 120061: iCiotti fc Ostrikerl l2007t ISiiacki et al l 
2007t iHopkins et al.ll2007bl )" which can result in the self- 



regulation of the gro wth of the black hole (see, e.g. 
IDi Matteo et al.l 120051 ). In this model we would expect 
black holes to grow rapidly during their early lifetime 
(i.e. while at low mass), until some point at which the 
black hole feedback begins to significantly affect its envi- 
ronment, resulting in a noticeable decline in growth rate. 
This effect has been observed in individu al black hole 
histo r ies, but such investigat ions (see, e.g. ISiiacki et all 
120091 : IDi Matteo et all 120111 ) have tended to focus on 
the largest mass black holes, primarily to explain how 
black holes could grow rapidly enough to produce the ex- 
tremely large masses (~ 10 9 M Q by z ~ 6) foun d in obser- 
vations by the Sloan Digital Sky Survey (e.g. iFan et all 
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20061 Uiang eTai]|2009] ) . In this paper we take advantage 
of a new, very large simulation to investigate the growth 
histories of early universe black holes across a wide range 
of masses, probing both the mean and the distribution of 
growth rates for black holes across a wide range of masses 
and luminosities, and provide fits for these distributions. 

2. METHOD 

In this paper we use a new cosmological hydrody- 
namic simulation of a 533 h^ 1 Mpc box specifically 
intended for high-redshift investigations. The simula- 
tion uses the massively parallel cosmolocial TreePM-SPH 
code P-GADGET (a n updated version of GADGET- 
2, see ISpringell 120051 ) incor porating a multi-phase ISM 
model with star formation (|Springel fc Hernquistl I2003T) 
and b l ack hole accretion and feedback (jSpringel et all 
I2005al IDi Matteo et all 120051 ). has a gravitational soft- 
ening length of 5 ft. -1 kpc and mass resolution of 2.8 x 
10 s M Q for dark matter and 5.7 x 1O 7 M for gas. 

Within the simulation, black holes are modeled as col- 
lisionless sink particles which form in newly emerging 
and resolved dark matter halos. These halos are found 
by calling a friends of friends group finder at regular in- 
tervals (in time intervals spaced by A log a = log 1.25). 
Any group above a threshold mass of 5 x lO lo ft- _1 Af 
not already containing a black hole is provided one by 
converting its densest particle to a sink particle with a 
seed mass of Mbh.scccI = 5 x 10 5 ft. _1 M Q . This seeding 
prescription is chosen to reasonably match the expected 
formation of supermassive black holes by gas directly col- 
lapsing to BHs with Mrh ~ Af seec j (e.g. iBromm fc Loebl 
120031 : iBegelman e~ al. 2006) or by Po pHI stars collaps- 
ing t o ~ W 2 M ( t ) BHs a t z ~ 30 (|Bromm fc Larsonl 
120041: lYoshida et ail 120061 ) followed by sufficient expo- 
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nential growth to reach M see d by the time the host halo 
reaches ~ 10 10 M Q . Following insertion, BHs grow in 
mass by accretion of surrounding gas and by merging 
with other black holes. Gas is accreted according to 

47rG 2 Mg H pBH 



M BH = 



( e 2 + ^2 )3 /2 ' , where p B H is the local gas density 
(determined from the gas particles within the black hole 
kernel), c s is the local sound speed, v is the velocity of 
the BH relative to the surrounding gas, and a is intro- 
duced to correct for the reduction of the gas density close 
to the BH due to our effective sub-resolution model for 
the ISM. To allow for the initial rapid BH growth nec- 
essary to produce sufficiently massive BHs at early time 
(~ 10 9 M Q by z ~ 6) we allow for mildly super-Eddington 
accretion, but limit it to a maximum of 3 x MEdd to pre- 
vent artificially high values. 

The BH is assumed to radiate with a bolometric lumi- 
nosity proportional to the accretion rate, L = t]Mbhc 2 
(jShakura fc Sunyaevl 119731 ). where the radiative effi- 
ciency 77 is fixed to 0.1 throughout the simulation and 
our analysis. To model the expected coupling between 
the liberated radiation and the surrounding gas, 5 per 
cent of the luminosity is isotropically deposited to the 
local black hole kernel as thermal energy. The 5 per cent 
value for the coupling factor is based on galaxy merger 
simulations such that t he normalization of th e Mbh — o" 
relation is reproduced (|Di Matteo et al.ll2005l ). 

The second mode of black hole growth is through merg- 
ers which occur when dark matter halos merge into a sin- 
gle halo, such that their black holes fall toward the center 
of the new halo, eventually merging with one another. In 
cosmological volumes, it is not possible to directly model 
the physics of the infalling BHs at the smallest scales, so a 
sub-resolution model is used. Since the mergers typically 
occur at the center of a galaxy (i.e. a gas-rich environ- 
ment), we assume the final coalescence will be rapid (e.g. 
IMaver et all 120071) , so we merge the BHs once they are 
within the spatial resolution of the simulation. However, 
to prevent merging of BHs which are rapidly passing one 
another, mergers are prevented if the BHs' velocity rela- 
tive to one another is too high (comparable to the local 
sound speed). 

The model used for black hole creation, accretion 
and feedback has been in vestigated and dis cusse d 



(2007); 



in ISiiacki et all 



|Pi Matteo et al 
Colberg fc di Matted "(120081); ISiiacki et al . ,_ 
DeGraf et all (|2010l) ; IDegraf et all (|2011l) . finding it 
does a good job reproducing the M bh — o relation 
the total black hole mass density (|Di Matteo et al.l 
l200l . the QLF (jPeGraf et a l.l l2010f). and the e xpected 
black hole clustering behavior iDegraf et al.l f20TTh . 
This simple model thus appears to model the growth, 
activity, and evolution of supermassive black holes in 
a cosmological context surprisingly well (though the 
detailed treatment of the accretion physics is infeasible 
for cosmological scale simu lations). We also note that 
iBooth fc Schavd (|2009t ) and I Johansson et al.l (j2008l ) have 
adopted a very similar model, and have independently 
in vestigated the paramet er space of the reference model 
of Di Mat teo et al.l (|2008f) , as well as varying some of the 
underlying prescriptions. In addition, this simulation 
has previously been used to inve stigate the growth of th e 
first very massive black holes dDi Matteo et al.l 1 201 lh , 
statistical properties of quasars ()DeGraf et al.1 120111) . 



and large scale high-resolution imaging (Fe ng et"al"1 
1 2 1 If) . For further details on the simulation methods 
an d convergence studies done for similar simulations, 
see lDi Matteo et al.l (|2008D . 

Because the simulation saves the complete set of black 
hole properties (mass, accretion rate, position, local gas 
density, sound speed, velocity, and BH velocity relative 
to local gas) for each BH at every timestep, the black hole 
output for such a large simulation is prohibatively diffi- 
cult to analyze using p revious techniques. For this rea- 
son, ILopez et al.l (|2011[ ) developed a relational database 
management system specifically for this simulation. A 
similar strategy has also b een followed in the analysis of 
the M illenium simulation (|Lemson fc Virgo Consortium! 
2006). In addition to providing a substantially more 
efficient query system for extracting information, this 
database is significantly more flexible than traditional 
approaches. For a complete sum mary of the databa se 
format and its efficiency, please see ILopez et al.l ([201 ID . 

3. RESULTS 
3.1. Typical Black Hole Growth Rates 
To quantify the growth rate of black holes, we use the 
mean Eddington ratio (A = ^ BH ) which we calculate 

-^edd 

for each black hole over a finite time interval. Because 
we have the complete BH growth history, we are able to 
compute this quantity based solely on the gas accretion, 
and neglect any mass gained through black hole mergers 
(though we find the mass gained by mergers to be small 
enough to have a negligible effect on our results) . In Fig- 
ureQ]we plot (A) /(l + z) 3 as a function of A/bh, initial for 
several redshift ranges. We plot (A) /(l + z) 3 rather than 
(A) for two reasons: First, to show that the dependence 
of (A) on Mbh is independent of redshift (at least for 
z > 4.75), and second to show that (A) oc (1 + z) 3 . 

Regardless of redshift considered, we find similar be- 
havior for Eddington ratio with respect to mass: more 
massive black holes grow faster than low mass black holes 
up to a peak growth rate at Mbh ~ 4 — 8 x 10 7 M Q , while 
the black holes above this characteristic mass grow more 
slowly. Thus black holes grow fastest (relative to their 
current mass) while at intermediate masses, and grow 
slower at higher mass. 

We find this peak in the Eddington ratio to be caused 
by the change in the local gas density available for fuel- 
ing BH growth. We plot the evolution in local gas den- 
sity (pbh, the density of gas contributing to Mbh) with 
mass in Figure [I] showing a clear peak at ~ 5 x 10 7 M Q . 
We note that neither the sound speed nor the BH veloc- 
ity (the other factors in the calculation of Mbh) exhibit 
a peak with respect to Mbh, confirming that the peak 
Eddington ratio is caused by the evolution in the local 
gas density. To show how the gas density evolves, in 
Figure [2] we show the gas density profiles around BHs 
below the Eddington ratio peak (~ lO 7 A/ - blue), at 
the peak (~ 5 x 10 7 Af Q - green), and above the peak 
(~ 4 x 10 8 M Q - red), each averaged across 100 BHs. In 
general we find the gas density profile to grow with ATbh 
until A/bh ~ 5 x 10 7 Af Q (as expected for BHs found in 
more massive halos). Above ~ 5 x 10 7 M Q the gas density 
away from the BH continues to grow, but the innermost 
density is suppressed, with the suppression growing with 
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Fig. 1. — Colored lines: The mean Eddington ratio ((A)) as a 



function of A/bh f° r several redshift ranges, scaled by 



(1+zY 



tdth 



Poisson error bars [Note that the datapoints' x-positions for each 
z-bin have been shifted to the right (3% increase for each z-bin) 
such that the error bars are distinguishable]. We also show the 
typical host halo mass corresponding to the given BH mass on the 
top axis. Filled circles: Average gas density at the BHs position 
(Pbh) f° r z = 4.75 — 5. 




r [kpc] 

Fig. 2. — Gas density profiles averaged among 100 black holes 
with mass ~ lO 7 Af (blue), ~ 5 X 1O 7 M (green), ~ 4 X 1O 8 M 
(red). Dotted line shows the gravitational softening length. 

Mbh in both magnitude and distance. This suppres- 
sion of the local gas density is caused by the feedback of 
the black hole, with the stronger feedb ack of high-mass 
BHs p roducing the strongest effect fsee lDi Matteo et"al"1 
()2011|) for detailed investigation of feedback among mas- 
sive BHs). 

We also show the typical mass of halos hosting a given 
Mbh on the top axis, noting that the Eddington ratio 
peaks at a host halo mass of ~ 7x 10 11 — lx 1O 12 M . This 
mass very closely matches the critical shock heating scale 
of ~ 6 x lO n M (|Dekel k Birnboiml [20061 IDekel et al.l 
120091 and consistent with our simulation), above which 
infalling gas is shock heated near the virial radius to the 



0.0 



-0.5 

A 

S 1 -1.0 



-1.5 
-2.0 



-*Shen et ol. 201 1 
MassiveBlack: 

M„„>10 7 M S T 
m ; < 21 * 




Pbh 



-22.0 



-22.5 ^ 
E 



-23.0 3 
a* 
o 

-23.5 



Fig. 3. — Redshift evolution of the Eddington ratio for black holes 
with Mbh > 1O 7 M (shaded region shows 1-cr standard deviation 
in log(A)) and i -band magnitude m,j < 21 (green line) compared 
with data from Shcn & Kcllv (2011) (black asterisks). We also 
show the evolution in the gas density around BHs for comparison 
(blue dashed line). 



virial temperature of the halo. IDekel k Birnboiml (|2006f ) 
suggest that in these halos AGN feedback becomes more 
significant, since the dilute shock- heated gas will be more 
susceptible to heating and pushing by the central AGN. 
This would thus produce a suppression in the gas density 
profile, consistent with the picture described above and 
the downturn in Figure [TJ 

In addition to the evolution in A with Mbh, Figure [1] 
also shows that A evolves with redshift as ~ (1 + z) 3 , 
which is also caused by the evolution in the local gas 
density. In Figure [3] we show the evolution in (log (A)) 
with redshift among A^bh > 10 7 M Q BHs (shaded region, 
showing l-cr standard deviation). We plot the average 
gas density at the BH (blue dashed line), showing the 
evolution in A is primarily caused by the evolution in pbh 
(recall Mbh oc p bh)- We also compar e to observational 
measurements of IS hen k Kellvl (|2011[) (black asterisks), 
showing that this general redshift evolution is consistent 
with current observations, and the normalization is ap- 
proximately consistent if we use a similar magnitude cut 
(i-band magnitude rrii < 21 - green line). 

3.2. Eddington Ratio Distributions 

In addition to investigating the mean Eddington ra- 
tio, we also study the distribution of A among com- 
parable BHs. Previous work on the A-distribution 
has often found rou ghly log-normal distributions using 
both observational (|Kollmeier et al.l [20061 INetzer et all 
[20071: INetzer fc Trakhtenbrotl 120071 : Iwillott et all 120101 : 
Trakhtcnbrot et al . 2011) and phen omenological ap - 
proaches (jShankar et al.ll2011l ) [though lAird et all (|2011f ) 
find A to follow a power law when selected for host stellar 
mass, rather than BH mass]. However, these observa- 
tional studies necessarily incorporate several uncertain- 
ties, such as sample selection and scatter in black hole 
mass estimators, which we can bypass, using our sim- 
ulation to probe our black holes' Eddington ratios di- 
rectly. In Figure 2] we show the distribution of Eddington 
ratios among black holes selected by Mbh (black his- 
tograms). We find that the distribution produced by 
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our simulation is ind eed log-normal, in keeping with ob- 
serva t ional findings (IKollmeier et aL 2006; N etzer et all 
20071 iNetzer fc Trakhtenbrotl l2007t Iwillott et al.l 120101: 
Trak htenbrot et al.ll2011h . In particular, we note that the 
distribution remains log-normal regardless of the mass 
considered, with Figure |4] showing this holds among black 
holes that are below, at, and above the peak observed in 
Figure [TJ 

Because we find A to follow a log-normal distribution 
and the mean of that distribution obeys a well-defined 
curve with Mbh (Figure [T]), we are able to provide a 
general fitting formula for P(A|Mbh, z), the probability 
distribution of black hole Eddington ratios as a function 
of redshift and black hole mass: 



10.00 



P(X\M BH ,z) 



(1) 



where fj, m and a m are the mean and standard deviation 
of ln(A), respectively, and are fit by a m ~ 0.39 and 



(2) 



with A ~ .00094, Mp = 5 x 1O 7 M , and cr - 0.85. In 
Figure |4] we plot the distribution predicted by Equations 
HB42] (red curve) compared the the actual distribution, 
showing that this simple formula is capable of reproduc- 
ing the distribution of A for BHs in our simulation across 
a wide range of masses and redshifts, without requiring 
knowledge of individual black hole environments. 

In addition to the distribution for a mass-selected sam- 
ple, in Figure [5] we show the Eddington ratio distribu- 
tion from our simulation (red histogram) comp a red to 
the observed distribution from IKollmeier et al.l (2006) 
(black histogram) for two luminosity selected samples. 
We again note that the distribution is described by a 
roughly log-normal distribution, and that our simulation 
is approximately consistent with observational results. 

Furthermore, by combining P(X\Mbh,z) with the 
black hole mass function ($bh) we can obtain the Ed- 
dington ratio probability distribution for a luminosity- 
selected sample: 



P(X\L BB ,z) 



$bh(M bh )P(\\M bh , z) 
Jo 00 $bh(M bh )P(A|Mbh, z)d\ 



(3) 



where A/bh = a 1 ^ n Figure [H we plot this pre- 

dicted probability distribution (using our simulation's 
mass function) in red, showing P(X\L BB , z) is well pre- 
dicted in this manner. We note that this approach is 
significant as it provides a potentially powerful tool for 
constraining the black hole mass function using observa- 
tions of the Eddington ratio distribution. We show this 
in Figure[5]by plotting P(A|L bh, z ) based on three differ- 
ent local mass functions: the lShankar et ail (I2009T) mass 
function (dashed green); t he iShankar et all ( | 2009) mass 
function derived from the iHopkins et al.l (|2007b|) lumi- 
nosity function (dash ed blue), and the mass function of 
IHopkins et al.l (|2007al ) (dashed pink) . Because P(A|L B h) 
is sensitive to the slope of 3>bh, the distribution of A at 
high Lbh (where the mass function is steepest) varies 
substantially with the mass function used, suggesting 
that with improved statistics from upcoming surveys, we 
could use the observed P(A|Lbh) to constrain the slope 




Fig. 4. — Eddington ratio distribution for black holes at three 
different mass scales (black) and the predicted distribution from 
Equations [T&rSl (red curves). 

of the black hole mass function at high redshift, even 
without measurements of the black hole masses. 

4. CONCLUSIONS 

With a new large-scale simulation, we show that the 
growth of black holes tends to follow a typical growth 
pattern. In general, we find that black holes grow more 
rapidly at higher redshift than comparable black holes at 
lower redshift, characterized by A oc (1 + z) 3 . This scal- 
ing is caused by the redshift evolution in the gas density 
about the black holes , and is comparable to current ob- 
servational data from [Shen & Kellv (201j]). 

The typical Eddington ratio also scales with Mbh such 
that A peaks at M B h ~4-8x 10 7 M q (typically found 
in halos of ~ 7 x 10 11 — 1 x 1O 12 M ). This peak is caused 
by evolution in the density of the gas at halo centers that 
is available to fuel black hole growth. In general, more 
massive black holes are found in more massive halos with 
correspondingly higher gas densities, hence A grows with 
Mbh for low masses. However, above Mbh ~ 5x 1O 7 M 
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black hole feedback has a sufficiently strong effect on the 
local environmen to suppress the density of the nearby 
gas. Thus although these more massive black holes are 
found in more massive halos with correspondingly higher 
gas densities in general, the feedback has significantly 
lessened the density of the innermost gas where accre- 
tion occurs. This suppression of the local gas density is 
stronger for more massive BHs, and causes A to decrease 
for M BH £ 5 x 10 7 M Q . 

Although the local environment is important for the 
accretion rate of individual black holes, we show that 
the distribution of Eddington ratios follows a roughly 
log- normal distribution regardless of the black hole pop- 
ulation considered, consistent with current observational 
findings. We use this, together with the evolution in (A), 
to provide a simple fitting formula for the distribution 
of Eddington ratio with (Mbh, z). This general forumla 
can be used for predicting the growth/evolution of black 
hole populations in theoretical and semi-analytic models 
(such as the evolution of the black hole mass function), 
for predicting the mass of observed high-redshift quasars, 
and, in conjunction with upcoming observations of the A- 
distribution, to constrain the slope of the high-redshift 
black hole mass function. 
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Fig. 5. — Distribution of Eddington ratios for BHs in our sim- 
ulation (red his togram) compared with observational data from 
IKollmeier et all 1120061 ) (black histogram) for two luminosity bins. 
We also show the predicted distribution based on our fitting func- 
tion _Q5qjrationsj3l using our simulation's mass function (solid red), 
the ISha nkar ct al. (200J3) base mass function (das hed green), the 
Shankar ct al. ( 2009) mass function derived from the Hopkins ct al. 
( 2007b ) luminosi t y funct ion (dashed blue), and the mass function 
of lHopkins et al.l |2007al ) (dashed pink). 
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